Effective continuous model for surface states and thin films of three-dimensional 

topological insulators 
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Two-dimensional effective continuous models are derived for the surface states and thin films of 
the three-dimensional topological insulator (3DTI). Starting from an effective model for 3DT1 based 
on the first principles calculation [Zhang et al, Nat. Phys. 5, 438 (2009)], we present solutions for 
both the surface states in a semi-infinite boundary condition and in thin film with finite thickness. 
The coupling between opposite topological surfaces and structure inversion asymmetry (S1A) give 
rise to gapped Dirac hyperbolas with Rashba-like splittings in energy spectrum. Besides, the SIA 
leads to asymmetric distributions of wavefunctions for the surface states along the film growth 
direction, making some branches in the energy spectra much harder than others to be probed by 
light. These features agree well with the recent angle-resolved photoemission spectra of Bi2Se3 films 
grown on SiC substrate [Zhang et al, arXiv: 0911.3706]. More importantly, using the parameters 
fitted by experimental data, the result indicates that the thin film Bi2Se3 lies in quantum spin Hall 
region based on the calculation of the Chern number and the Z2 invariant. In addition, strong SIA 
always intends to destroy the quantum spin Hall state. 

PACS numbers: 



I. INTRODUCTION 

Topological insulators (TIs), which are band insula- 
tors with topologically protected edge or surface states, 
have attracted increasing attention recently 1 . A well- 
known paradigm of topological insulator is the quantum 
Hall effect, in which the cyclotron motion of electrons 
in a strong magnetic held gives rise to insulating bulk 
states but one-way conducting states propagating along 
edges of system 2 .. The idea was generalized to a graphene 
model with spin-orbit coupling, which exhibits the quan- 
tum spin Hall (QSH) state 3,4 . Later, the realization 
of an existing QSH matter was predicted theoretically 5 
and soon confirmed experimentally^ in two-dimensional 
(2D) HgTe/CdTe quantum wells. Furthermore it was 
found that the QSH state can be induced even by the dis- 
orders or impurities 8 - - — . Meanwhile, the concept was also 
generalized for three-dimensional (3D) TIs, which are 
3D band insulators surrounded by 2D conducting surface 
states with quantum spin texture 1 ^—. Bi^Sbi-^, an al- 
loy with complex structure of surface states, was hrst con- 
firmed as a 3DT I 15 i 16 . Soon after that it was verified by 
both experiment a 17 i 18 and first-principles calculations 19 
that stoichiometric crystals Bi2X3 (X=Se, Te) are TIs 
with well-defined single Dirac cone of surface states and 
extra large band gaps comparable with room tempera- 
ture. The Dirac fermions in the surface states of 3DTI 
obey the 2+1 Dirac equations and reveal a lot of uncon- 
ventional properties and possible applications, such as 
the topological magneto-electric effect^ and Majorana 
fermions for fault-tolerant quantum computing 2 ^—. 

Thanks to the state-of-art semiconductor technologies, 
low-dimensional structures of Bi2X3 can be routinely fab- 
ricated into ultra-thin film o 27 ' 28 and nanoribbons 2 ^. This 
stimulates several theoretical works on the thin films of 
3DTIs— ~— . For further studies of the transport and op- 



tical properties of 3DTI films and their potential appli- 
cations in spintronics and quantum information, it is de- 
sirable to establish an effective continuous model for thin 
films of TIs. 

In this paper, we present an effective continuous model 
for the surface states and ultra-thin film of TIs. Start- 
ing with a 3D effective low-energy model based on the 
first principles calculations 1 ^, we first present the solu- 
tions for the surface states and the corresponding spectra 
for a semi-infinite boundary condition of gapless Dirac 
Fermions and for the thin film of TIs. The finite size 
effect of spatial confinement in a thin film leads to a 
massive Dirac model which may exhibit the QSH ef- 
fect. Within the same theoretical framework, a structure 
inversion asymmetry (SIA) term is further introduced 
in this work to account for the influence of substrate, 
providing a description of the Rashba-like energy spec- 
tra observed in the angle resolved photoemission spec- 
tra (ARPES) in the recent experiment on Bi2Se3 films^ 8 -. 
We derived the parameter conditions for the formation 
of QSH effect in a thin film in the absence and presence 
of the SIA. By analyzing the fitting parameters with the 
help of the Chern number and the Z2 invariant, we iden- 
tified the ultrathin films of Bi2Se3 in the QSH phase in 
the experiment. 

The paper is organized as follows. In Sec. |TT]we intro- 
duce an anisotropic 3D Hamiltonian for 3DTI, which is 
a starting point of the present work. With this Hamilto- 
nian, we present detailed solutions to the thin film in two 
different boundary conditions. In Sec. IIII1 effective con- 
tinuous models are established for the surface states and 
thin film of 3DTI. Within the framework of this effective 
continuous model, the structure inversion asymmetry is 
taken into account and an effective Hamiltonian for SIA 
is derived in Sec. IIV1 In Sec. [V] we apply the model to 
newly fabricated thin film Bi2Se3 and demonstrate that 
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thin films of Bi2Se3 are in the QSH regime. Finally, a 
conclusion is presented in Sec. IVII 



II. MODEL AND GENERAL SOLUTIONS FOR 
3DTI 

A. Model for 3DTI 




FIG. 1: Schematic of a topological insulator film grown on 
substrate. The grown direction is defined as z axis. The 
thickness of the film is L. 

As shown in Fig. [TJ we will consider a thin film grown 
along z direction. The thickness of the film is L. We 
assume translational symmetry in x-y plane so that the 
wave numbers k x and k y are good quantum numbers. We 
start with the effective model proposed to describe the 
bulk states near the T point for the bulk Bi2Se3^. The 
states are mainly contributed by four hybridized states 
of Se and Bi p z orbitals, denoted as {|pl+,t), |f>2j,t}) 
|pl+,|), \p2~,l)}, where + (— ) stands for the even (odd) 
parity. The Hamiltonian is given by 



H{k) = e (k)7 4 



M(k) 

-iAxd z 


A 2 k+ 



-M{k) 
A 2 k+ 






A 2 k- 
MQt) 
iA x d z 



k x ±ik y , e (k) = C-Dxd 2 



A 2 k- 


iA x d z 
-M(k) 

(1) 

D 2 k 2 ,M(k) = 



where k± 

M + Bidl-B 2 k 2 , and k 2 = k%+k%, with A 1 , A 2 , B 1 , B 2 , 
C, Di, D 2 , and M the model parameters. This model has 
the time reversal symmetry and the inversion symmetry. 
Though we start with a concrete model, the conclusion 
in this paper should be applicable to other topological 
insulator films. We shall demonstrate that this model 
for the bulk states can produce the surface states with 
appropriate boundary condition. 



B. General solutions of the surface states 



putting a four-component trial solution 

ijj = ijj x e Xz 



(2) 



into the Schrddinger equation [E is the eigenvalue of en- 
ergy) 



H(k,-id z )xP = Ei>, 
the secular equation 

det\H(k,-i\)-E\ = 



(3) 



(4) 



gives four solutions of X(E), denoted as (3\ a (E), with 
a G {1,2}, G {+,-}, and 



X a (E) 



F 



+ (-1) 



a-l 



R 



2D+D_ v ' 2L>+LL 
where for convenience we have defined 



(5) 



A\ +D+(E- L x ) + D_(E- L 2 ), 



F = A 2 

R = F 2 - AD + D- [(E — Li)(E - L 2 ) — A^k+k-], 

D± = Di±B u 

In = C + M+(D 2 -B 2 )k 2 , 

L 2 = C -M+{D 2 + B 2 )k 2 . 



(6) 



Because of double degeneracy, each of the four (3\ a (E) 
corresponds to two linearly independent four-component 
vectors, found as 



tpa/31 



i>afi 2 = 



D+X 2 a -L 2 + E 
-iAi{p\ a ) 


A 2 k+ 

A 2 k_ 


iA x {p\ a ) 
D-Xl-Li+E 



(7) 



(8) 



The general solution should be a linear combination of 
these eight functions 



$(E,k,z) 



^EE C a ^ aM e^ z , (9) 



ct=l,2 p=± 7=1,2 



with the superposition coefficients C Q ^ 7 to be determined 
by boundary conditions. In the following, we will con- 
sider two different boundary conditions: one is semi- 
infinite focusing on only one surface at z = - 0; the other 
includes two opposite surfaces at z = ±L/2. In both 
cases we assume open boundary conditions (\P = 0) for 
the surface states at the surfaces. 



Following the method by Zhou et ali^, the general 
solution for either the bulk states or the surface states 
can be derived analytically. Despite the existence of time- 
reversal symmetry, the A 2 k± term couples opposite spins 
in Hamiltonian (1), and one has to solve a 4 x 4 matrix, 
instead of the simplified 2x2 one in the 2D case^. By 



C. Solutions for the surface states with 
semi-infinite boundary conditions 

The surface states have a finite distribution near the 
boundary. For a film thick enough that the states at op- 
posite surfaces barely couple to each other, we can focus 
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on just one surface. Without loss of generality, we study 
a system from z = to +00. The boundary condition is 
given as 



*(z = 0) = and *(z -> +00) = 0. 



(10) 



The condition of ^l{z — > +00) = requires that \P con- 
tains only the four terms in which (3 = — and the real 
part of A a is positive. 

Applying the boundary conditions of Eq. (fT0|) to the 
general solution of Eq. ((9]), the secular equation of the 
nontrivial solution to the coefficients C a ^ 7 leads to 



(Ai+A 2 ) 5 



(11) 



which along with Eq. ([5]) gives the dispersion of the sur- 
face states 



E± = C 



DiM 

By 



(ft - 



(12) 

Near the T point, the dispersion shows a massless 
Dirac cone in k space, with the Fermi velocity vf — 

{A 2 /K)^l - (§^) 2 , instead of plain A 2 /H as in RefJ£. 

The wave functions for E± are found as 



m_ = c Q _ 




-A+z 



-a; 



(13) 



where A„ are short for X a (E = E±) according to Eq. 
(|5|), taxiip = ky/k x , and C± are the normalization factors. 
The properties of the solution to A Q determine the spatial 
distribution of the wave functions. Generally speaking, 
the edge states exist if Ai and A2 are both real or com- 
plex conjugate partners. For either case, there should be 
inequality relations 



M 



> 0, 



Bi 

D+D- < 



(14) 



The edge states distribute mostly near the surface (z = 
0), with the scale of the decay length about Aj~ 2 for real 
A12 or [Re(Ai, 2 )] _1 for complex Xx 2 . In the former case, 
the wavefunctions decay exponentially and monotonously 
away from the surface (not from z = 0); while in latter 



case the decaying is accompanied by a periodical oscilla- 
tion, which can be easily seen from the wavefunctions in 
Eq. (HH). In addition, there exist complex solutions to 
A a when 



-D+D- 



< 



AM 



(15) 



D. 



Solutions for finite-thickness boundary 
conditions 



When the thickness of the film is comparable with the 
characteristic length 1/Ai^ of the surface states, there is 
coupling between the states on opposite surfaces. One 
has to consider the boundary conditions at both surfaces 
simultaneously. Without loss of generality, we will con- 
sider the top surface is located at z = L/2 and the bottom 
surface at —L/2. The boundary conditions are given as 



$(z = ±-) = 0. 



(16) 



In this case, the general solution consists of all eight 
linearly independent functions. Applying the boundary 
conditions in Eq. (fTB"]) to the general solution of Eq. © , 
the secular equation of the nontrivial solution to the su- 
perposition coefficients C a /3j leads to a transcendental 
equation 



tanh 



A 2 L 

2 



tanh 



2 



A1A2 



tanh 



2 



tanh 



\ 2 L ■ 
(17) 

In a large L limit, tanh reduces to 1 , then Eq. (fTT)) 
can recover the result in Eq. ([IT]). With the help of Eq. 
©, Eq. (IT7|) can be used to identify the energy spectra 
and the values of A Q numerically. 

Due to the finite size effect, the coupling between 
the states at the top and bottom surfaces will open an 
energy gap22r— . We define the gap as A = E + — E- at 
the T point, where E + and EL are two solutions of Eq. 
(fTT]) . For X a L > 1 and A 2 > Ai (L can be finite), the 
approximate expression for A can be found. If X a is real, 
the gap can be approximated by 



A\A 1 D+D_M\ 



-Ai L 



(18) 



which decays exponentially as a function of L. Fig. HJa) 
shows the gap as a function of thickness, in which a set of 
model parameters used to fit the ARPES of 4QL Bi 2 Sc 3 
thin film is employed, as listed in the first row of Tab. UJ 
For some other materials there may exist complex Ai = 
A 2 and we can define Ai = a — ib and A 2 = a + ib, where 
a > 0, b > according to Eq. ([5]). In this case, the gap 
is found as 



A 



8\AiD+D-M\ 



^J-BHA^+AD+D.M) 



e- aL sin(bL), (19) 
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III. 



EFFECTIVE CONTINUOUS MODELS 




FIG. 2: (Color online) [(a)(b)] The energy gap A = E+ - E- 
and [(c) (d)] the model parameter B [defined in Eq. (|42[l ] as 
functions of the film thickness L. Circles correspond to the 
numerical results of the transcendental equations (|30|l and 
(I31|l . Solid and dash lines correspond to the approximate 
formulas to A when L is finite [Eqs. (|18l) and (|19|) ] or very 
small [Eq. (|5U|l ]. respectively. All the parameters are adopted 
(a) by fitting experimental results of 4QL Bi2Se3, and (b) 
from the numerical fitting for the first principles calculation 
of Bi 2 Se^, as listed in Tab. HI 



TABLE I: Two sets of parameters for the 3D Dirac model. 
The first row is extracted from our effective model parameters 
for 4QL Bi2Se3 film in table[II] and the second row is adopted 



from the first principles calculation 



M A x A 2 Si 


B 2 


C 


Di 


D 2 


(eV) (eVA) (eVA) (eVA 2 ) 


(eVA 2 ) 


(eV) 


(eVA 2 ) 


(eVA 2 ) 


0.28 3.3 4.1 1.5 


-54.1 


-0.0068 


1.2 


-30.1 


0.28 2.2 4.1 10 


56.6 


-0.0068 


1.3 


19.6 



The solutions of the surface states and thin film of 
3DTI can be applied to calculate physical properties ex- 
plicitly. For instance, we can see whether the ground 
state of a thin film exhibits QSHE or not by calculating 
the Chern number or Z 2 invariant. It is also desirable 
to establish an effective continuous model to explore the 
properties of these surface states especially when other 
interactions have to be taken into account. For this pur- 
pose, in this section we derive an effective low-energy and 
continuous models for the surface states and thin film of 
3DTI. 

Due to the low-energy long-wavelength nature of the 
Dirac cone of the surface electrons, we can use the solu- 
tions of the surface states at the T point as a basis to ex- 
pand the Hamiltonian H(k) in Eq.(l), which will be valid 
when the energy is limited within the band gap between 
the conduction and valence bands. This is equivalent to 
a truncation approximation as we exclude the solutions 
for the bulk states in the basis. In this approach, the 
Hamiltonian in Eq. ([!]) can be expressed as 



H(k) = H Q (k = 0)+AH, 



where 

with 
h(Ai, 

and 



h{A x ) 
/i(-Ai) 



(22) 



(23) 



AH = 



-D_d 2 z +C + M -iAid z 
-iAid z -D+dl + C-M 



D 2 k 2 - B 2 k 2 a z A 2 k-a x 
A 2 k + cr x D 2 k 2 - B 2 k 2 a z 



(24) 



(25) 



The first term can be solved exactly, and the last term 
describes the behaviors of electrons near the T point. 



with 



At 



A. Basis states at V point 



2 ^-D+D_ ' 




(20) 
(21) 



According to thi s result , the oscillation period of the gap 
ir/b becomes tt^JBi/M when A\ — 0, in accordance with 
the result obtained by Liu et a£22. Fig. [U[b) shows the 
gap oscillation by using the model parameters listed in 
the second entry of Tab. |TJ Besides, the sine function 
implies that A may be negative. Later we will see that 
the sign of A can be found by solving E+ and E°_ from 
Eqs. (|3"tj|) and (|3"Tj) . respectively. 



Hq in Eq. (|23j) is block-diagonal. Its solution can be 
found by solving each block separately, i.e., h(Ai) 1 ^^ = 
£;^ t and h(-Ai)^^ = E^^. Because the lower block 
is the "time" reversal of the upper block, the solutions 
satisfy ^4.(2) = Q^f(z), where 9 = —ia y JC is the 
time-reversal operator, with a y the y component of the 
Pauli matrices and K, the complex conjugation operation. 
Equivalently, we can replace A\ by — A\ in all the results 
for the upper block, to obtain those for the lower block. 
Therefore, we only need to solve h(Ai). Following the 
same approach in Sec. HH we put a two-component trial 
solution 



1> 



t 



(26) 
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into 



h(A 1 ,-id z )^ = E^, 



(27) 



the secular equation for a nontrivial solution yields four 
roots of X(E), denoted as (3\ a , with (3 e { + ,— } and 
a G {1,2}. Note that here X a is short for A Q (fc = 0) 
in Eq. ([5]). Each f3X a corresponds to a two-component 
vector 



D+\ 2 a -l 2 + E 



(28) 



The general solution is a linear combination of the four 
linearly independent two-component vectors 



E E Cartttef**: 

a=l,2 13=+,- 



(29) 



Applying the boundary conditions Eq. (|T6|) to this gen- 
eral solution, we obtain two transcendental equations, 



and 



(C-M-E- L>+A 2 )A 2 _ tanh(^) 
(C - M - E - D+X^Xj. ~ tanh(^)' 



(C-M-E- D + Xl)X 2 _ tanh(^) 
(C - M - E - D + X 2 2 )X! ~ tanh(^)' 



(30) 



(31) 



The solutions to Eqs. (|30l) and (j3"Tj) give two energies at 
the T point, designated as E+ = E + (k = 0) and E°_ = 
E^(k = 0), respectively. The eigen wavefunctions for E+ 
and E°_ are, respectively, 



t 



x(Ai) = = a 



iAxft 

-D+ri2f+ 
iAxfZ 



(32) 
(33) 



cosh(Aiz) cosh(A 2 z) 



where C± are the normalization factors. The superscripts 
of f± and r]f 2 stand for E± , and the subscripts of f± for 
parity, respectively. The expressions for f± and r]f 2 are 
given by 

/+(*) = 
f-(z) = 



cosh(^) C osh(^; 
sinh(Aiz) sinh(A 2 ;z) 



sinh(^) sinh(^) 
(Ai) 2 - (A 2 ) 2 



Ai coth(^) - A 2 coth(^p) 
(Ai) 2 - (A 2 ) 2 



E=E ( ' 



Ai tanh(^) - A 2 tanh(^) 



(34) 
(35) 
(36) 
(37) 



The energy spectra and wavefunctions of the lower block 
h(— Ai) of H can be obtained directly by replacing A-y 
by —A\. Based on the above discussions, the four eigen- 
states of Hq can be given by 



<p(Ai) 


, * 2 = 


x(Ai) 













<f(-Ax) 



, $4 = 







X(-At) 



(38) 



with $i — > $3 and $ 2 — > &4 under the time-reversal op- 
eration. We should emphasize that these four solutions 
are for the surface states, and the solutions for the bulk 
states are not presented here. We use the four states as 
the basis states, and other states are discarded (except 
that in Fig. [3] where four extra bulk states are also in- 
cluded by the same approach), because of a large gap 
between the valence and conduction bands. 



B. Effective model for 3DTI films 

With the help of the four states Eq. ([35} at the V 
point, we can expand Hamiltonian Eq. (jTJ) to obtain a 
new effective Hamiltonian 



f L/2 

H cS = / dz[$ 1 ,$ 4 ,$ 2 ,$3] t ^[$l,$4,$ 2 ,$ 3 ], (39) 
J-L/2 



where for convenience, we organize the sequence of the 
basis states following {$i, $4, 4> 2 , $3}. Under the reor- 
ganized basis, the effective Hamiltonian is found as 



h+ 
h- 



(40) 



with 



h+ = E — Dk 
h_ = E - Dk 2 



A 

2 



Bk 2 



A* 2 k+ - 

-%+Bk 2 
-A 2 k+ 



A 2 k_ 
f+Bk 2 



^ 

±-Bk 2 



(41) 



and 



B = 



Bi — B 2 



D = 



B\ + B 2 



2 2 

E = (E+ + E_)/2, A = E+ 

Bi = Ba{(p(Ai)\<r M \<p(Ai)), 
B 2 = B 2 ( x (A 1 )\a z \x(Ai)), 
A 2 = A 2 ( V {A 1 )\a x \ X {-A 1 )). 



-D 2 , 
-E_, 



(42) 



We find that A 2 here can be either real or purely imagi- 
nary (see Appendix I VIII for details), classifying the model 
into two cases: 
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Case I is for a real A 2 = Hvp, and the effective Hamil- 
tonian is further written as 

h Tz — E — Dk 2 + Hvft z (t -k + T z a z ( — — Bk 2 ), (easel); 

(43) 

and case II for a purely imaginary A 2 = ihv-p, 

h Tz = E — Dk 2 + hvF((7xk) z +T z a z (^-~Bk 2 ), (case II) 

(44) 

with t z — ±1 corresponding to the upper (lower) 2x2 
block in Eq. (|40|) . vf the defined Fermi velocity and a 
and k here only refer to the components in x-y plane. In 
fact, these two effective Hamiltonian can consist of the 
invariants of the irreducible representation D 1 / 2 of SU(2) 
groups. 

Eq. ([4"Tj) can also be expressed in terms of the Pauli 
matrices 

h Tz =E - Dk 2 + d ■ a, (45) 



with 



T z (hi>Fk x , hvpky, — — Bk 2 ), (easel) 
(HvFky, —hvFk x , t z ( — — Bk 2 )). (case II) 

(46) 

d(A:) vectors in case I and case II, respectively, corre- 
spond to Dresselhaus- and Rashba-like textures. Note 
that case I is essentially the effective 4x4 model for the 
CdTe/HgTe quantum wells^. However, we find that case 
I only occurs for quite a small range of thickness. For 
most thicknesses of interest, A 2 is pure imaginary. There- 
fore, we only focus on case II in the following discussions. 
By far, we have reduced the anisotropic 3D Dirac model 
into a generalized effective model for 2D thin films, under 
the freestanding open boundary conditions. 



which has the same dispersion as Eq. (fl"2j) and the 

same Fermi velocity vf — 4f yl ~~ (§t) 2 as f° r tne 
semi-infinite boundary conditions. In an isotropic case, 
D\ = D2 and B\ = B2, the quadratic term disappears 
and we have a linear dispersion for the Dirac cone. Fi- 
nally it is noticed that the models for the surface states at 
the top and bottom surface have the same form assumed 
\ a L ^> 1. We will see that these results work well even 
for films down to five quintuple layers (QL) of atoms in 
thickness (1 QL is about 1 nm). 



D. The ultra-thin limit 

Another opposite limit is L — > 0, which is a little bit 
complicated since X a L does not approach to zero when 
L is very small. In Eq. (13TH) . the left side has an order 
of L 2 when L — > 0, so tanh(^p) must have the order of 
L~ 2 , which means 



tanh(-^) 



=s> Ai = i-. 

Li 



(49) 



Combining this result with Eq. ([5]), the model becomes 



Dm 2 



D 2 k 2 +A 2 (ax k) z +T z ( 



B l7 r 2 



B 2 k 2 )a x 



(50) 

It is found that a finite energy gap opens at k = 0, i.e., 
A = 2Bitt 2 /L 2 as shown in Fig. [2] Note that this result 
in the L — > limit even provides a rough estimate of the 
gap for most thicknesses. Besides, the continuum limit 
generally assumed in this work also works well even for 
only several quintuple layers. 



IV. STRUCTURE INVERSION ASYMMETRY 



A. Structure Inversion Asymmetry 



C. Effective continuous model for surface states 

Despite the simple explicit form, the parameters in 
Hamiltonian (|40[) need to be determined numerically. Be- 
fore that, we can take two limits to see their behaviors. 
The first limit is X a L 3> I, for a = 1,2. In this case, 
tanh(^) ±2 1, and both Eqs. flU and reduce to 

(C-M-E-D+XDXz = (C-M-E-D+Xj)Xi. (47) 

Solving this equation, we have an effective continuous 
model for the surface states (ss) of 3D topological insu- 
lator as 



C 



DxM 
Bi 



{D 2 -B 2 ^)k 2 



-A 2X /l-(^na x k y -a y k x ), (48) 



A recent experiment revealed that the substrate on 
which the film is grown influences dramatically electronic 
structure inside the film. Because the top surface of the 
film is usually exposed to the vacuum and the bottom 
surface is attached to a substrate, the inversion symmetry 
does not hold along z direction, leading to the Rashba- 
like energy spectra for the gapped surface states. In this 
case, an extra term that describes the structure inversion 
asymmetry (SIA) needs to be taken into account in the 
effective model. 

We use the same method as that in Sec. IIIII to in- 
clude the SIA term. Without loss of generality, we add 
a potential energy Viz) into the Hamiltonian. Generally 
speaking, V{z) can be expressed as V(z) = V s (z) + V a (z), 
in which V s (z) = V s (-z) and V a (z) = -V a (-z). The 
symmetric term V s could contribute to the mass term A 
in the effective model, which may lead to an energy split- 
ting of the Dirac cone at the T point. We do not discuss 
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it in details in this paper. Here we focus on the case 
of the antisymmetric term, V{z) — V a (z), which breaks 
the top-bottom inversion symmetry in the Hamiltonian. 
A detailed analysis demonstrates that V a {z) couples $1 
($3) to $2 ($4), which can be readily seen according 
to their spin and parity natures. The modified effective 
Hamiltonian in the presence of V{z) becomes 



H, 



SIA 



0^0 
V* 
V * 
V 



(51) 



where 



r= / L/2 dzMA^VaizMA,)). 

l-L/2 



(52) 



Comparing this definition with that of A\ in Eq. (|42p. we 
find that V also can be either real or purely imaginary. In 
the case of a purely imaginary (case II) A2, V must be 
real (see Appendix IVII[) , and the effective Hamiltonian 
with SIA can be written as 



rrSIA 



h+(k) 

Va 



Va 
h-(k) 



(53) 



In the case of a real A2 , V must be purely imaginary, and 
the effective Hamiltonian with SIA then has the form 



rrSIA 



h+(k) 

-Va z 



Va z 
h_{k) 



(54) 



Without the SIA term, the effective Hamiltonian (|44| 
gives the energy spectra of the gapped surface states as 



E ± =E Q - Dk 2 ± 



(— -Bk 2 ) 2 
y 2 ' 



(hv F ) 2 k 2 , (55) 



where + (— ) sign stands for the conduction (valence) 
band, each of which has double spin degeneracy due to 
time-reversal symmetry. When the SIA term is included, 
the Hamiltonian (ISTj) gives 



given with the help of Fig. [4] On the left is for a thicker 
freestanding symmetric TI film, and it has a single gap- 
less Dirac cone on each of its two surfaces, with the solid 
and dash lines for the top and bottom surface, respec- 
tively. The two Dirac cones are degenerate. The top of 
Fig. @] indicates that the inter-surface coupling across 
an ultrathin film will turn the Dirac cones into gapped 
Dirac hyperbolas. On the bottom of Fig. 01 SIA lifts 
the Dirac cone at the top surface while lowers the Dirac 
cone at the bottom surface. The potential difference at 
the top and bottom surfaces removes the degeneracy of 
the Dirac cones. On the right of Fig. @J the coexistence 
of both the inter-surface coupling and SIA leads to two 
gapped Dirac hyperbolas which also split in fc-direction, 
as shown in Fig. [3J 



Top surface 




FIG. 3: (Color online) Energy spectra of surface states (four 
branches in the middle) and several branches of bulk states 
(those at the top and the bottom) for a film with the thickness 
of 4 QL in the presence of the structure inversion asymmetry. 
The color of lines corresponds to the spatial distribution of 
the wavefunctions in z direction. Dark blue (light green) rep- 
resents that the wavefunctions mainly distribute on the side of 
the top (substrate) surface. The model parameters are listed 
in the first row of Table HI 



E h ± = E -Dk 2 ± J(^-Bk 2 )2 + (\V\ + hv F k)i, 



E 2 ,± = E -Dk 2 ± ] J(j-Bk 2 )2 + (\V\-hv F k)2, 

(56) 

where the extra index 1 (2) stands for the inner (outer) 
branches of the conduction or valence bands. The en- 
ergy spectra in the presence of V is shown in Fig. [3J 
Each spin-degenerate dispersion in Eq. (l55|) shifts away 
from each other along k axis. Both the conduction and 
valence bands show Rashba-like splitting. An intuitive 
understanding of the energy spectra in Fig. [3J can be 



B. Location of the surface states 

Location of the surface states can be revealed by eval- 
uating the expectation of position z of these states. The 
spatial distributions along the z direction of a state ip a 
can be evaluated by the expectation of position in z di- 
rection (z), 



Z\%p a \ 2 dz. 



By this definition, (z) a e [— |r] 
for a symmetric spatial distribution 



and 



(57) 



becomes 
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inter-surface 
coupling 

A ■ 



Dirac cones 




SIA 
+ 

inter-surface 
coupling 



SIA 



FIG. 4: (Color online) The evolution of the doubly-degenerate 
gapless Dirac cones for the 2D surface states, in the pres- 
ence of both the inter-surface coupling and structure inver- 
sion asymmetry, into the gapped hyperbolas that also split 
in fc-direction. The blue solid and green dashed lines corre- 
spond to the states residing near the top and bottom surfaces, 
respectively. 



With the SIA, the eigen wavefunctions are found as 



i(t + E l £)k 
+ hv F k)k+ 
i(\V\ +hv F k)k 



(58) 



with Ef = E 1± -E + Dk 2 , t= f - Bk 2 , and, 



2± 



2 y /E% lt (E^ t +t)k 



(\V\ - hv F k)k + 
i(-\V\ + hv F k)k 
(t 



E^ ut )k A 



(59) 



with £^ ut = E 2 ±-E + Dk 2 . Fig. [3] demonstrates (z) by 
the brightness of lines, with dark blue for (z) = -| (the 
top surface), and light green for (z) = — ^ (the substrate 
or bottom surface). 

For a thin film of 4 quintuple layers (QL), L = 3.8nm, 
it is found that the two surface states are well separated 
and dominantly distributed near the two surfaces. The 
averaged (z) ~ which is about 2/3 of a QL (w L/6) 
deviating from the surface. In this case the top and bot- 
tom surface states are well defined even without the SIA 
(V = 0). The average value remains almost unchanged 
in a large range of k. However, the crossing point of the 
spectra of the top and bottom surface states, the aver- 
aged (z) changes from +L/3 to 0, and then goes to the 
value of —L/3. This demonstrates that the finite thick- 
ness makes the two states couple with each other as their 
wave functions along the z direction have a finite overlap. 



As a result the two states open an energy gap as in the 
case of edge states in QSH system^. The value of the 
gap is a function of L as shown in Fig. 2(a) and (b). 
Near this region, (z) varies from (z) ~ L/3 to —L/3, and 
becomes zero exactly when two states are mixed com- 
pletely. For a large L, we find that the averaged distance 
of the surface states deviating from the surface remains 
about 1 QL. 

Simply speaking, the states close to the top surface 
are easier to be probed by light than those close to the 
bottom surface. This provides a hint to understand why 
there are branches in energy spectra with much more 
faint ARPES signals^,. 



V. THIN FILM Bi 2 Sc 3 AND QSH STATES 

In this section, we will investigate the realization of 
QSH effect in thin films and apply the effective model to 
the thin film Bi2Se3. When the system does not break the 
inversion symmetry, the effective Hamiltonian is block- 
diagonalized by t z — ±1. This is in a good agreement 
with the theory by Murakami et al2£. In this case we 
can define a t z -dependent Chern number (Hall conduc- 
tance) for each block like the spin Chern number—, from 
which the nontrivial QSH phase can be identified. Af- 
ter introducing the SIA term, the r z -dependent Chern 
number loses its meaning as the two blocks are mixed to- 
gether. However, we can still employ the Z2 topological 
classification 4 , which requires no inversion symmetry, to 
identify possible QSH thin films in experiment. 



A. QSH effect without SIA 

Considering the block-diagonal form of the effective 
model without SIA (|40|) . we can derive the Hall conduc- 
tance for each block, separately. For the 2x2 Hamilto- 
nian in terms of the d(fc) vectors and Pauli matrices in 
Eq. (|45|) , the Kubo formula for the Hall conductance can 
be generally expressed a o 37 ' 38 



( fk- - fk,+) dd a ddp 



> xy 



2m v 

k 



where is the volume of the system, d the norm of 
{d x ,d y ,d z ), / fc ,± = l/{cxp[(£ ± (fc) - ti)/k B T] + 1} the 
Fermi distribution function of electron (+) and hole (— ) 
bands, with \x the chemical potential, ks the Boltzmann 
constant, and T the temperature. 

At zero temperature and when the chemical potential 
[i lies between the band gap (—3^, 3 )' ^ ne Fermi func- 
tions reduce to fk,+ — and fk t - — 1. In this case we 
have 3 ^ 



-r»2^[sga(A) + sgn(B)]. 



(61) 



This result intuitively shows that only when B and A 
have the same sign, the Chern number is equal to +1 or 
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— 1, which is topologically nontrivial, and the Hall con- 
ductance is quantized to be ±e 2 /h. In other words, the 
QSH depends not only on the sign of A at the T point 
but also on that of B for k large enough. Experimen- 
tally, the ^-dependent Hall conductance can be probed 
by the nonlocal measurement, just like that for the 2D 
CdTe/HgTe quantum wells^. 



B. QSH effect with SIA: Z 2 invariant 



(a) V = 



(b) V = 0.053 eV 




FIG. 5: (Color online) [(a) and (b)] P(k) for four combina- 
tions of A and B, in absence (a) and the presence (b) of a 
small SIA V. |A|=0.07 eV and |B| = 10 eV A 2 , (c) The con- 
tour C is used to count the number of the pairs of the zeros 
of P(k), which form a circle ring when AB > 0. 



with i,j run over all the bands below the Fermi surface, 
i.e., "01- and ip2- m the present case according to Eqs. 
([H} and ([59]) . Based on the spin nature of the basis 
states {$i, $4, $2, $3} in our effective model, the time- 
reversal operator here is defined as = ia x ®a y K, where 
<7 X and <j y are the x- and y-component of Pauli matrices, 
respectively, and K. the complex conjugate operator. The 
number of pairs of zeros can be counted by evaluating 
the winding of the phase of P(k) around a contour C 
enclosing half of the complex plane of k = k x + ik y , 



iS\. 



(64) 



Because the model is isotropic, we can choose C to enclose 
the upper half plane, the integral then reduces to only 
the path along k x -axis while the part of the half-circle 
integral vanishes for 6 > and |k| — > +00. 

In the absence of the SIA term, P(k) is found for the 
Hamiltonian (|4"4")) as 



P(k) = 



-%+Bk 2 



(A _ Bfc 2)2 + (^ p )2 fc 2 



(65) 



in which one can check that the zero points exist only 
when k 2 = A/2B > 0, and form a circle ring. Along 
fcz-axis only one of a pair of zeros in the ring is enclosed 
in the contour C, which gives a Z2 index 1 = 1. This 
defines the nontrivial QSH phase, and is in consistence 
with the conclusion by the Hall conductance in Eq. f|0 1 [) . 

In the presence of a small SIA term V < hv F a/|A/2B], 
with the help of the eigen wavefunctions ({58} and (|59[) . 
real P(k) can be found (after a U(l) rotation) as 



P{k) = 



(t + E™)(t + E°_ ut ) + \V\ 2 - (hv F k) 
2^E in E_ out (t + E™){t + E° ut ) 

f sgn(Hv F k - |% A > 
I 1, A<0 ' 



(66) 



In the presence of SIA, V couples the blocks h+ and 
so the ^-dependent Hall conductance becomes non- 
sense. Following Kane and Mele^, we can employ the Z2 
topological classification to give a criterion of the QSH 
phase, because it does not require inversion symmetry as 
a necessary condition. The Z2 index can be obtained 
by counting the number of pairs of complex zeros of 
P(k) = PfL4(k)], where the Pfaffian is defined as 



Pf[A(k)] 



1 



2 n n\ 



E 



Pcrmutationsof {i± . . . . .i2 n } 



■At 



(62) 



where the sgn is to secure the continuity of P(k). One 
can check that P(0) = — sgn(A) and P(oo) = sgn(S) 
. Besides, for a small V, the behavior of P(k) between 
P(0) and P(oo) will not change qualitatively (see Fig. 
|£J}. Therefore, for AB > 0, P(k x ,0) should still have odd 
pairs of zeros. For a large V > Hvf\/\A/2B\, 



P(k) = 



(t + E™)(t + E° ut ) + \V\ 2 - (hv F k) 2 
2 x /E™E° ut (t + E™)(t + E° ut ) 



x f sgn(hv F k - \V\), B < 
I I, B>0 



(67) 



in which N counts the number of times of permutations, 
and A(k) is a 2n order anti-symmetric matrix defined by 
the overlaps of time reversal 



A^(k) = ( Ui (k)|0| Uj (k)) 



(63) 



One can check for this case P(0)P(oo) is always positive 
thus P(k) has even pairs of zeros, regardless of the signs 
and values of A and B. In other words, a large SIA will 
always destroy the quantum spin Hall phase. 



C. Thin film Bi 2 Se 3 and QSH effect 



TABLE II: Fitting parameters to the Bi2Se3 thin films, using 
the energy spectra Eq. (156[l from our effective model, [adopted 
fromRef.-^] 



Layers 


E 


D 


A 


B 


Wp 


\v\ 


(QL) 


(eV) 


(eVA 2 ) 


(eV) 


(eVA 2 ) 


(10 5 m/s) 


(eV) 


2 


-0.470 


-14.4 


0.252 


21.8 


4.47 





3 


-0.407 


-9.7 


0.138 


18.0 


4.58 


0.038 


4 


-0.363 


-8.0 


0.070 


10.0 


4.25 


0.053 


5 


-0.345 


-15.3 


0.041 


5.0 


4.30 


0.057 


6 


-0.324 


-13.0 








4.28 


0.068 



Recently, thickness-dependent band structure of 
molecular beam epitaxy grown ultrathin films Bi2Sc3 
was investigated by in-situ angle-resolved photoemission 
spectroscopy^ 8 -. An energy gap was first observed exper- 
imentally in the surface states of Bi2Se3 below the thick- 
ness of six quintuple layers, which confirms theoretical 
prediction as a finite size effector—. 

Table [TT] shows the fitting parameters to the ARPES 
data of Bi 2 Se3 thin films 28 using the energy spectra for- 
mula [Eq. (|56j) ]. For the films with thickness ranging 
from 2 QL to 5 QL, all of them satisfy sgn(A£>) > and 
V < hv-py / \A/2B\, hence the films are possibly in the 
QSH regime. We identify that only 2 QL, 3 QL, and 4 
QL belong to the nontrivial case for potential QSH effect. 
5 QL is an exceptional case that the fitted parameters B 
and D do not satisfy the existence condition of an edge 
states solution^. The condition of B 2 < D 2 will lead to 
the band gap closing for a large k. However, it is under- 
stood that the model is only valid near the T point, and 
the fitting parameters are limited to the case of small k. 
And the band gap was measured clearly for the film of 
5QL. 

It was previously predicted, using the parameters from 
the first-principles calculation^, that the gap A should 
oscillate as a function of the film thickness^— . However, 
this oscillation is not reflected in the measured results. 



D. QSH effect of SI A and the edge states 

In the quantum Hall effect the Chern number of the 
bulk states has an explicit correspondence to the num- 
ber of edge states in an open boundary condition^. In 
topological insulator or QSH system, the Z 2 topological 
invariant has also a relation to the number of helical edge 
states^. As a supplementary support to the above con- 
clusion, we demonstrated the presence of edge states in 
a periodic boundary condition along the x-direction and 
an open boundary condition (say along y-direction) im- 
posed in a geometry of strip of the thin film by means 
of numerical calculation. Using the parameters in Table 



10 

I, we have concluded that a stripe of 2 - 4 QL will ex- 
hibit helical edge states. More specifically, we present the 
energy dispersion for 4QL in Fig. [5] There is a doubly- 
degenerate Dirac point inside the gap of the 2D surface 
states for 4 QL in consistence with the results obtained 
in the above sections. 



4QL 




MA" 1 ) 

FIG. 6: (Color online) The energy spectra by the tight bind- 
ing calculation for 4 QL of Bi2Se3 thin films with the width 
along y-direction 200 nm. The parameters are given in table 



VI. CONCLUSIONS 

We derived two-dimensional effective continuous mod- 
els for the surface states and thin films of three- 
dimensional topological insulators (3DTI). A gapless 
Dirac cone was confirmed for the surface states of a 3DTI. 
For a thin film, the coupling between opposite topologi- 
cal surface states in space opens an energy gap, and the 
Dirac cone evolves into a gapped Dirac hyperbola. The 
thin film may break the top-bottom symmetry. For ex- 
ample, the thin film grows on a substrate, and possesses 
the structural inversion asymmetry (SIA) . This SIA leads 
to a Rashba-like coupling and energy splitting in the mo- 
mentum space. It also leads to asymmetric distributions 
of states along film growth direction. 

The ARPES measurements on Bi2Sc3 films have 
demonstrated that the surface spectra opens a visible en- 
ergy gap when the thickness is below 6QLsj2£ The energy 
gap was observed to be a function of the thickness of thin 
film, and in a good consistence with theoretical prediction 
as a finite size effect of the thickness of thin film. The 
Rashba-like splitting was measured clearly in the thin 
film of 2 to 6 QLs. This can be explained very well from 
the inclusion of the SIA. Since the thin film was grown 
on a SiC substrate and the other surface is exposed to 
the vacuum this fact results in the SIA in the thin film. 
Another direct evidence to support the SIA is the sig- 
nal intensity pattern of the energy spectra of ARPES. 
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Usually the surface states are located dominantly near 
the top and bottom surfaces. The signal intensity for 
these two branches of energy spectra of ARPES are dif- 
ferent. The SI A will cause the coupling between two 
surface states near their crossing point. That is why the 
Rashba-like splitting of the ARPES spectra has a bright 
crossing point near the T point, with one branch bright 
and the other almost invisible. Thus the SIA term can be 
used to describe the ARPES measurements on the thin 
film Bi 2 Se3 very well. 

Our effective model demonstrates that the 3DTI can 
be reduced to an two-dimensional quantum spin Hall sys- 
tem due to the spatial confinement. Strictly speaking, 
the system is no longer a 3DTI in the original sense once 
the energy gap opens in the surface bands, since the Z2 
invariant for the bulk states becomes zero. However the 
surface bands themselves may contribute a non-trivial 
one in the Z2 invariant even when the SIA term is in- 
cluded. Our calculation demonstrates that a strong SIA 
always intends to destroy the quantum spin Hall effect. 
A critical value for SIA exists, at the point there is a 
transition from a topological trivial to non-trivial phases. 
Based on the model parameters fitted from the experi- 
mental data of ARPES, we conclude that the thin film 
Bi2Se3 should exhibit quantum spin Hall effect once the 
energy gap opens in the surface spectra due to the spatial 
confinement of the thin film. 
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In this appendix, we demonstrate that both the pa- 
rameters Ai and V in the effective model (|40| can be 
either real or purely imaginary, and the product of A 2 V 
must be pure imaginary. By putting the wavefunctions 



TABLE IV: Four possible groups of E+ and E- , and resulting 
values of A2 and V . 



£4 






E- 




A 2 


V 


case 


A 




case 


A 


imaginary 


real 


case 


A 




case 


B 


real 


imaginary 


case 


B 




case 


A 


real 


imaginary 


case 


B 




case 


B 


imaginary 


real 


Eqs. 


©2 


| and 


(ESI) 


into the definitions 


in Eqs. (02) and 


(El 


we 


have 











A 2 = iAiAzD+C+C- J L dz{ri2ft*f++vf*f+*f-} 

V = C + C_ [ 2 dzV{z)[D\4*^ fl* 7+ + A 2 J+* fZ]. 

(68) 
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VII. APPENDIX: MODEL PARAMETERS A 2 
AND V 



TABLE III: Four possible combinations of Ai and A2 , accord- 
ing to Eq. ([5]), and resulting f± and 771,2 according to Eq. 
(|38[1 . According to Eq. ©, A^ < A2, so there does not exist 
a case when Ai is real and A2 is pure imaginary. 





Ai 


A 2 


f± m,2 


case A 


A 2 


A! 


imaginary real 




real 


real 


real real 


case B 


imaginary 


real 


real real 




imaginary 


imaginary 


real real 



For arbitrary energy, Eq. ([5]) requires that the values of 
Ai and A2 can only be one of the combinations shown in 
Tab. III. By putting these combinations into Eq. (|M)l - 
P7p. one can show that for the first entry, 771 and 772 are 
real while /+ and /_ are pure imaginary, referred as the 
case A; while for entries 2-4, all of 771, 772, /+, and /_ 
are real, referred as the case B. For an arbitrary group of 
E+ and E_ , each of them belongs to either the case A or 
B, leading to four possibilities, as shown in Tab. IIVI In 
particular, according to Eq. (fT5| . when A\j{— D+D-) > 
4M/B\, there is no complex Ai^, corresponding to the 
last row of Tab. IIVI i.e., A2 is pure imaginary while V is 
real. 
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